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Abstract 

Many biological and physical systems exhibit population-density dependent transitions to syn- 
chronized oscillations in a process often termed "dynamical quorum sensing". Synchronization 
frequently arises through chemical communication via signaling molecules distributed through an 
external media. We study a simple theoretical model for dynamical quorum sensing: a heteroge- 
nous population of limit-cycle oscillators diffusively coupled through a common media. We show 
that this model exhibits a rich phase diagram with four qualitatively distinct mechanisms fueling 
population-dependent transitions to global oscillations, including a new type of transition we term 
"dynamic death". We derive a single pair of analytic equations that allows us to calculate all 
phase boundaries as a function of population density and show that the model reproduces many 
of the qualitative features of recent experiments of BZ catalytic particles as well as synthetically 
engineered bacteria. 



Unicellular organisms often undertake complex collective behaviors in response to envi- 
ronmental and population cues. A beautiful example of this phenomenon is the population- 
density dependent transition to synchronized oscillations observed in communicating cell 
populations recently termed dynamical quorum sensing pp. Density-dependent synchro- 
nization has been observed in a wide variety of biological systems including suspensions of 
yeast in nutrient solutions [2], starving cellular colonies of the social amoeba Dictyostelium 
[3], and synthetically engineered bacteria [I]. Such transitions have also been observed in 
experimental studies of electrochemical oscillators and Belousov-Zhabotinsky (BZ) catalytic 
particles [5J E] . 

Previous theoretical work has shown that oscillators coupled through quorum sensing 
can display synchronized oscillations [THS]. Recently, a dynamic quorum sensing transi- 
tion was found [2] in a simple model of coupled identical limit-cycle oscillators introduced 
to study synchronization in yeast populations. Additionally, numerical studies of BZ cat- 
alytic particles indicate that heterogeneity in oscillator populations leads to interesting new 
phenomenon [HI [TO] - The study of population-density dependent synchronization in het- 
erogeneous populations of oscillators is still in its infancy, in stark contrast to oscillators 
with direct coupling where many analytic results are available [T0l - fl2"] . 

In this paper, we consider a large population of limit-cycle oscillators with a distribution 
of natural frequencies, coupled diffusively through a common external media. Our work 
generalizes earlier models [2] and exhibits extremely rich dynamics as the coupling strength, 
population density, and frequency distribution are varied. We derive several analytic results 
and find that model exhibits a rich phase diagram. In particular, there exist four quali- 
tatively different mechanisms leading to synchronization as a function of coupling strength 
and population density: (1) a Kuramoto-like incoherence to coherence transition, (2) an 
amplitude death transition due to oscillator heterogeneity, (3) a loss of both global and in- 
dividual oscillations due degradation in the external media, and (4) a new type of transition 
we term "dynamic death" where the external media dynamics are not fast enough to support 
global oscillations. We show that the model reproduces many qualitative features observed 
in recent experiments on heterogeneous populations of BZ catalytic particles [3] as well as 
synthetically engineered bacteria jlj. 

To illustrate this diverse set of phenomena, we introduce a simple model of iV 3> 1 coupled 
limit-cycle oscillators where the amplitude and phase of individual oscillators are represented 
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by a complex number Zj, (j = 1...N), with natural frequency Uj. The oscillators are 
diffusively coupled to an external media, represented by a complex number Z, through a 
coupling D. Furthermore, when chemicals leave the oscillators and enter the medium, they 
are diluted by a factor a = Vi n t/V ex t <C 1, which is the ratio of the volume of the entire 
system to the that of an individual oscillator. The external media is degraded at a rate J. 
The dynamics of the system are captured by the equation 

dz- 

= (A + iujj - Izjl^Zj - D{ Zj - Z) 

f t =aD^~Z)-JZ 

j 

where the Uj are drawn from a distribution h(oj) which we assume to be an even function 
about a mean frequency ojq. By introducing a dimensionless density, p = aN, and shifting 
to a frame rotating with frequency Uq, we can rewrite the equations above as 

dz- 

= (A + iujj - \zj\ 2 )zj - D(zj - Z) 
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where the frequencies Uj are now drawn from an even distribution g(oo) with mean zero. 

To build intuition for the system, it is helpful to consider the special case of a homoge- 
nous population where g(u) is a delta function and all u) 3 ; = in (JT1). This model was used 
previously [2] to model dynamical quorum sensing in yeast suspensions. For homogenous 
populations, the equations for all the Zj are identical and there are two possible behav- 
iors. The individual oscillators are quiescent with Z = Zj = or there are synchronized 
oscillations. We can compute the stability of the oscillator death state by linearizing the 
system around Zj = Z = and computing the eigenvalues, /i, of the corresponding lin- 
earized system. In this calculation, since all of the oscillators are identical, the dynamics 
are completely specified by two differential equations, one for the mean-field parameter 
z = Zj = Zj and one for Z. Oscillator death is stable when Re(/i) < for all eigen- 

values. The corresponding requirement that the trace be negative implies D > Xq in the 
oscillator death phase. Furthermore, the characteristic equation for the eigenvalues takes 
the form (// + + B + iu ) = pD 2 , with B = Dp + J and A = D - A . To find the 
phase boundary, we plug in p = a + ib and separate the characteristic equation into real and 



FIG. 1. (Top) Phase boundaries for homogeneous oscillators in the p vs. D plane for various 
values of J, ojq = 1, and Ao = 1. (Bottom) Heat map of the amplitude of collective oscillations 
from simulations of N = 40 identical oscillators, showing the transition from a "dynamic" death 
phase with J = to synchronized oscillations for D > 1. Inset: Real part of z and Z during low 
density oscillations showing ~ 1000-fold slowing of oscillations (D=2.5, p = 0.001). 



imaginary parts, 



a + B= ;°\: A l (2) 
(a + A) 2 + b 2 w 

(a + A) 2 + b 2 



This allows us to solve for b as a function of a, b(a) = — 2a+B+A an< ^ P^ U S ^ ms ^ n ^° (22). The 



resulting equation can be analyzed graphically plotting the left and right sides of (22) as a 
function of a [13]. Since the characteristic equation is quadratic, there are two solutions, a 
solution with negative a which guarantees the stability of the external medium, and a second 
solution which can change sign depending on parameters. Examining the aforementioned 



plot, it is clear that if the left-hand side of (22) is greater than the right-hand side at a = 



then the second solution must also be negative. Thus, the amplitude death phase is stable 
when 

(Dp + J)(D-X ) (pD + J + D-\ ) 2 

P D 2 ~ {pD + J + D- A ) 2 + tu 2, U 

where we have rewritten A and B in terms of the original parameters of the model. 

This equation indicates that there are two qualitative ways to stabilize amplitude death. 

First, when J ^> 1, the left side is much larger than the right, indicating that oscillations 

are lost due to degradation of the external medium. Second, even when J = 0, if the 



natural frequency of the oscillators is large, uq 3> 1 or the density of oscillators is small, 
amplitude death can be stable because the right hand side of ^ is extremely small. The 
underlying reason for this is that since D—X > 0, isolated individual oscillators are silent and 
synchronization can only occur by transmitting information through the external medium. 
The external medium has an effective time scale, (p_D) _1 , on which it can respond. If uq 
is large, the medium cannot track the dynamics of the individual oscillators and the death 
phase is stable. We term this "dynamic death" to indicate that the underlying cause for 
the loss of oscillations are the dynamical properties of the external medium. Figure [3] shows 
the phase boundaries in this model as a function of J, p, and D. Notice the dynamic death 
region for small J and p. We have also confirmed the existence of this phase using numerical 
simulations (see Figure [3j bottom). 

We now analyze (jl]) for the case where the natural frequencies tOj are drawn from an even 
distribution g(u) with zero mean. In this case, the system has three phases: an amplitude 
death phase where all oscillators are quiet; global, synchronized oscillations; and an incoher- 
ent phase where individual elements are oscillating but the oscillations are unsynchronized. 
For finite densities, we did not numerically observe any unsteady behavior analogous to that 
seen in [10J. The stability boundary of the amplitude death phase can again be calculated 
as in the homogenous case by linearizing the system around the death state Zj = 0, Z = 
(for details see [13J). In the large N limit, the boundaries of stability for the death phase 
are D — A > and 

pD + J f D - A 



J dug{uJ \D-\ Q y + (b-uy (5) 

du) SMrFi x , 2 , 77 , x 2 - (6) 



pD 2 J " v J (D- X^ + {b-ujf 

It is useful to consider various limits of these equations. Notice that when g(cu) = S(cu), these 



equations reduce to (22) with a = as expected. Alternatively, consider the case p — > 00. 



In this limit, the left hand side of (24) is zero implying that 6 = 0, since g(u>) in an even 



function. Plugging this into (|5| yields a single equation for stability of the death state, 

pD + J f D-Xp 

duj 9{^)TT, x N 2 , ,. 2 - ( 7 ) 



This equation was derived in [TOl [TT] for the stability boundary of the death phase in a 
system of directly coupled limit-cycle oscillators. This follows naturally by noting that in the 
limit p — > 00, the external media can respond infinitely quickly. Thus, Z ext is equal to the 



order parameter of the system, Z ext = ^2j z j> an d the model reduces to the one studied 
in [TUHT2]. In this limit, the loss of oscillations is due to the heterogeneity of individual 
oscillator frequencies and has been termed amplitude death. These two limits show that a 
single set of equations (j5])-(24), capture all four qualitatively distinct types of transitions to 



synchronized oscillations. 

To gain further insight into the system, it is useful to consider the mean-field equations 
for the system. To do so, we put Zj = rje ldj and Z = Re 1 ^ into Q and equate real and 
imaginary parts: 

(IT ■ 

-± = (A - D - rf)rj + DR cos (0 - 9j) 



d6i DR 



— — = ujj H sin (0 — 6j. 

at ta 



% = ^ + f±'i*n«-S,). (8) 
We look for uniform rotating solutions whose angular frequency in the lab frame is Uq + b 



by requiring ^ = = and -£ = = b in (|32|). In general, these equations cannot be 



solved analytically. However for the special case R = corresponding to the death phase 



one finds that the equations reduce to (24) (see [13J). This allows us to attach a physical 
meaning to the parameter b. When oscillators are directly coupled to each other (p — > oo), 
they rotate at the mean frequency ojq and 6 = 0. In contrast, when oscillators are coupled to 
each other through the external media, there is an effective "viscosity " which slows down 
the oscillations so that they rotate with an angular frequency u + b, with b < 0. Thus, b 
measures the change in angular frequency of oscillations due to delays induced by the ex- 
ternal medium (see Figure [3] inset). Furthermore, increasing J decreases the absolute value 
of b and hence increases the angular frequency. Thus, somewhat surprisingly, the system 
exhibits positive period-amplitude coupling. This behavior was observed in a population of 
synthetically engineered bacteria in recent experiments though interestingly, the same 
phenomena occurs as degradation is decreased. The effect of time delays on synchroniza- 
tion of directly coupled oscillators was studied previously and the equations governing the 
stability of amplitude death bear some similarity to those found in this work [HI [15] . 
When D < X Q , the system can be incoherent, where individual oscillators are oscillating in 



an unsynchronized fashion. The stability equations for the incoherent phase were calculated 



by generalizing the calculations in [TU]. We looked for solutions of (32) of the form R = 0, 
r| = Ao — D, and 9j = u>jt. For such solutions, individual oscillators oscillate at their natural 
frequencies but there are no coherent oscillations. We calculated the stability boundary for 
incoherence by checking the stability of the state to small perturbations (see [13]). We 
find that for distributions g(oj), there is a tri-critical point on the line D = Ao where the 
incoherent phase, the synchronized oscillation phase, and the death phase meet. 

The derivation presented above is for arbitrary g(u). When g(u) is either a rectangular 



or Lorentzian distribution, we can perform the integrations in (24) explicitly (see [13] )• The 
answers are particularly simple for a Lorentzian distribution, g(cu) = ~ r^+oj 2 > anc ^ are gi ven 
by 

pD + J _ D-\ + T 
pD 2 ~ (D-X + T) 2 + b 2 

h + luJ * = b (o) 

pD 2 (D - A + r) 2 + b 2 ' y } 

These equations are identical to the homogenous case, (|22j), except with A — > A — T on the 
right hand side. Thus, heterogeneity "renormalizes" the distance individual units are from 
their Hopf bifurcation. An analogous set of equations-albeit more unwieldy-can also be 
derived for a rectangular frequency distribution, and Figure [2] shows the phase boundaries 
for this case as a function of p, D, and J. As expected, for D > Ao, the death phase 
and synchronized oscillations are both possible. For large D, as density is increased across 
the transition, the amplitude of the synchronized oscillations rises sharply with density. 
For smaller D, this rise in amplitude is less pronounced. When D < X , one also sees a 
Kuramoto-like transition from incoherent to synchronized oscillations. The same qualitative 
behavior was observed in recent experiments on BZ catalytic particles with a distribution 
of natural frequencies [3 |6] . 

In this letter, we considered the physics of limit-cycle oscillators diffusively coupled 
through an external media. We have found that there are four qualitatively different types of 
density-dependent transitions to synchronized oscillations including a new type of transition 
we term "dynamic death" , where the dynamics of the external medium are too slow to sup- 
port oscillations. This simple model captures many qualitative features seen in a variety of 
experiments on oscillators coupled diffusively through an external medium. For example, it 
was previously argued that when all oscillators are identical, the model is a good description 
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FIG. 2. (Top) Phase boundaries in p vs D plane for heterogeneous oscillators drawn from a 
rectangular distribution ( V = 0.6, ujq = 1, and A = 1 and J = .05 — 0.5 from bottom to top). 
(Bottom) Heat map of the amplitude of collective oscillation from simulations of N = 100 oscillators 
drawn from a rectangular frequency distribution with parameters above for J = 0.5. Notice the 
transitions from an incoherent phase with D < 1 to synchronized oscillations or death. Incoherence 
phase (spotted red region for D < 1) can be identified by the 1/iV 1 / 2 fluctuations of the order 
parameter. 



of glycolitic oscillations in suspensions of yeast cells. The model also shows how large am- 
plitude oscillations can emerge as one varies the density and how this behavior crosses-over 
into a Kuramoto-like transition as D is decreased (see Fig. [2|. These qualitative features 
are in good agreement with [HI E]- Finally, the model also captures many of the mean-field 
properties of coupled synthetically-engineered bacteria, including the sudden emergence of 
oscillations and scaling of amplitude and period of oscillations as one changes the external 
degradation rate J. However, unlike [I], the period and amplitude of the oscillations de- 
crease with increasing J. This discrepancy likely arises from the highly non-linear nature 
of the "degrade-and-fire" oscillations [TB] characterizing the synthetic bacteria. Despite this 
discrepancy, our results suggest that properly constructed simple models may be able to 
capture interesting, qualitative behaviors of coupled oscillators. In the future, it will be 
interesting to directly relate this simple model to more detailed models [TB] and extend 
the simple mean-field model of dynamical quorum sensing explored here to include spatial 
effects. 
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I. SUPPORTING INFORMATION: EQUATIONS FOR STABILITY OF AMPLI- 
TUDE DEATH 



In this section, we derive equations for the stability of the amplitude death phase. We 
start with the equations 



dz- 

= (A + iujj - \zj\ 2 )zj - D( Zj - Z) 
^ 1 ^ y^X z 3 - Z) - (J + iu )Z. 



dt 



N 



(10) 
(11) 



where we have transformed to a frame rotating at the mean frequency, ujq. To calculate the 



stability of amplitude death, we linearize these equations around the solution Zj 



Doing so yields the equations, 

d5zj 

~dT 
d5Z 

~df 



(Ao + iujj — D)5zj — D5Z 



These equations can be written in matrix form as 

( Sz x \ ( (Xo-D + icoi) • ■ 

(A — D + iuj 2 ) ■ ■ 



5zo 



5z N 
\SZ J 












-D 

-D 



= 0. 

(12) 
(13) 

\ fSzA 

Sz 2 



pD 

N 



pD 

N 



(Xo-D + ico N ) -D 

-{pD + J + iu ) J 



pD 
N 



5z N 
\5Z J 



(14) 



Denote the matrix on the right hand side of (14) by M for notational convenience. Stability 



requires the eigenvalues, /i, of M to satisfy Re[fi\ < 0. First notice that stability requires 
Re[Tr(M)] < 0. This gives the condition 



(A„-Z»-^>0. 



(15) 



We can also calculate the eigenvalues using the characteristic equation of the matrix, 
Det(fil — M) = 0. A straightforward calculation yields 

pD z 



(p+( P D + J + iu )) -(X -D + 

3=1 



IQJj) — 



e n (^-( a c 

s=l j=l,jj=s 



D + iujj 



(16) 
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In order to take the thermodynamic limit, we rewrite this equation as 



12 N 



( , + (pD + J + ^ )) = ^__±__ (17) 



In the thermodynamic limit N — > oo but with p held fixed, (15) and (17) become, respec- 
tively, 

(Xo-D)<0 (18) 

and 

P + {pD + J + iu )) = P D 2 [ du ( , 9{Uj) n ^ . r , (19) 

J p - (A - D + iu) 

where we have replaced the sum by an integral over the distribution function g(u) for the 
oscillator frequencies. In practice, it is often helpful to write this as two real equations. 
Substituting p = a + ib yields the coupled equations. 

a + pD + J f a + D — \ 



, dwq(oo)- -—7; — (20) 

pD 2 J yK ' (a + D - A ) 2 + (6 - u) 2 V ; 

" / duJ 9(u)—— ; jr ; (21) 



pD 2 J av (a + D — A ) 2 + (6 — a;) 2 
Stability requires that all solutions of these equations obey a < 0. By considering the mean- 
field equations derived below (and considering various limiting cases of the equations above), 
it is clear that the stability boundary can be found by putting a = in the above equations, 
i.e. there exists at most one solution with positive real part. 

We illustrate this explicitly for the case of uniform frequencies, where all the oscillators 
are identical g{oS) = S(u) (since by assumption the center of the distribution is around 
u> = 0). Then, the equations above can be rewritten as 

a + B^ ^^l (22) 
(a + A) 2 + b 2 v ' 

b + tu = . ~ P A D 9 b - . (23) 
(a + A) 2 + b 2 v ' 

where A = D — Xq and B = Dp + J. Dividing the equations and solving for b in terms 
of a to get b(a) = — 2a+B+A resu lt s i n an equation for a that can be analyzed graphically. 



Plotting the left-hand and right-hand sides of equation (23) in Figure |3| we see that there 
always exists a solution with negative real part. Recall that the characteristic equation is 
quadratic, assuring at most two solutions. The second solution can be either positive or 
negative and, since there are no other solutions with positive real part, setting a = safely 
identifies the phase boundary. 
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FIG. 3. Graphical analysis of the structure of the solutions to eq. (22) and (23). There always 
exists a solution with a < 0, shown as a black dot. A second solution, shown as a red dot, can 
have positive or negative a, and thus solving for a = properly identifies the phase boundary. 



Putting a = for the case with heterogeneous frequencies results in a pair of coupled 
integral equations that determine the boundary of stability of the death phase: 



P D + J f A I \ D ~ X ° (OA\ 



pD 2 J " v ' (D - A ) 2 + (b - uf 

du >g{u)T^ x , 77 — (25) 



pD 2 7 av (£) - A ) 2 + (b — lo) 2 



We can consider various limiting cases of the equations above. For uniform frequencies, 
the resulting equations are precisely those derived by Del Monte and co-workers. We can also 
consider the p — > oo limit. This corresponds to the case where the oscillators are directly 
coupled through the mean-field parameter. In this case, the left-hand side of the second 



equation in (25) is zero. Since g(oj) is an even function, this implies 6 = 0. Plugging this 



into the top equation in (25) and taking the limit p — > oo yields the equation 



h = h^ \ D -^ - (26) 



This is precisely the equation derived by Mirollo et al. for direct all-to-all coupling. 
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II. MEAN-FIELD EQUATIONS FOR FREQUENCY LOCKING 



Here we derive the mean-field equations for the model and show that they reproduce the 
results from the determinant calculations. We begin again with the equations 

dz- 

= ( Ao + iu _ \ z .\*) z . _ D ( Zj _ Z ) (27) 

f = ^E^'- Z )-( J +^) Z - ( 28 ) 

In writing these equations we have assumed that Uj are drawn from a normalized, even, 
probability distribution g(u>j). Now define Zj = rje l9j and Z = Re 1 ^. Then, we can rewrite 
these equations in polar coordinates to get 

dv ■ 

-2 = (A - D - rj) rj + DR cos (0 - 6j) (29) 

d6j DR , , „ . . 

-£ = Ui + — tintt-Oj) (30) 



dR _ pD 
~dt ~ ~W 



N 



^2 r j cos (0 - 9 j)-(pD + J) R (31) 



3=1 

£ = + (32) 

3=1 



We now look for uniform, rotating, locked solutions where ^ = — and ^ = = 6. 
In this case, the position of each oscillator is determined purely by its frequency, so we can 



regard each oscillator as a function of u. Using equations (29) and (30) and plugging in the 
expressions for the derivatives on the left hand side one can easily show that 

((w - b) cot (9 - 0) + A - D)(l + cot 2 (0 - <j>)){u - b) 2 = D 2 R 2 (33) 

Now, in the thermodynamic limit plugging in the desired solutions yields 

R(pD + J)=pD J duog{uj)r{uj) cos (0 - 6(u)) (34) 

b + u = pD J dug{u) r -^-sin((j)-6(u)), (35) 

where we have written r{oS) and 9(ou) to emphasize the amplitude and phase of each oscillator 

ing equation 

DRsm(9(u)-^) 
(u — b) 
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is a function of only the frequency. Using equation (30) and the fact that -nr = b yields 



Plugging this into (34) and (J35J) gives the equations 

sin (9(lj) — 0) cos (0(cj) — 0) 



(pD + J) 
pD 2 
b + u 
pD 2 



dug(oo) 



dug{ui 



sin (6{u) - 0) sin (0(w) - 0) 



(37) 
(38) 



Combined (37), (38), and (33) define the mean-field equations for the system for frequency 
locking with amplitude R. 



In general, we cannot solve this analytically because (33) is a cubic equation in cot. 



However, for the special case R = (oscillator death), we have the unique solution to (33) 
that 

tan(0(u;)-0) = ^^-. (39) 



An — -D 



Plugging this into the equations (34) and (35) yields the equations for amplitude death to 
be stable 



(pD + J) 
pD 2 

(b + u ) 
pD 2 



dujgioj) 



D-Xr 



dug(u 



(D - A ) 2 + (b - uj) 2 
(b-u) 



(40) 
(41) 



(D - An) 2 + (b - OJ) 2 

These are precisely our equations obtained from the derivative expansion. It gives as an 
interpretation of the imaginary part of the eigenvalue b as the shift of the frequency of the 
collective oscillations from uq. This should be contrasted with the case where oscillators are 
directly coupled and b = 0. 

Finally, we can ask intersect the stability boundary D = An. We can do this by taking 
the limit (D — A) — > in the equations above. A straightforward calculation shows that the 
equations reduce to 

(pD + J) 



pD 2 
(b + cup) 
pD 2 

where the P denotes the principal value. 



Tigib) 
P 



du 



9{u) 
uj — b 



(42) 



III. MEAN FIELD EQUATIONS FOR INCOHERENT STATE 



In this section, we derive the mean field equations for the stability of the incoherent state 
when D < A. In this case, we look for solutions of (29)-(32) of the form R = 0, r 2 = An — D. 
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In such a solution, individual oscillators oscillate at their natural frequencies but there is no 
coherent oscillations. 

We now examine the stability of the coherent state. To do so, we follow the analysis 
developed by Mathews et al. Define a density function p(r, 9, u, t) so that the fraction of 
oscillators of frequency u between r and r + dr and between 9 and 9 + d6 is prdOdr. The 
evoloution for p is given by the continuity equation 

dp 



Of 



+ V • (pu) = 



(43) 



where V is the velocity of oscillators given by v = (r,r#). Substituting (29) and (30) gives 



Q0 1 Q 1 (J 

H + -— (p [r 2 (a 2 -r 2 ) + K Rr cos(6 - (p)]) + -— (p [rcu - KRsm(9 - 0)]) = 0, (44) 



dt r dr 

where a 2 = (Ao — D). In the incoherent state 

5(r — a) 



P 



(45) 



We now consider a small perturbation in the radial and angular directions and check 
when the density is stable to these perturbations. In particular, consider 



p = 5(r-a-er 1 (6,u J t)) — + ef 1 {9,u,t) . 



1 



2nr 



(46) 



For such a perturbation, by the chain rule we have 



dr\ dn 



(47) 



Writing R = eR±, substituting in (|29|) and (|30|), and keeping terms first order in e yields 

(48) 



- 2aVi + DRt cos(# - 



dri dr\ 



89 dt' 

We seek solutions in which R\ and r\ are proportional to e^ x+tbS)t and we find that r\ must 
obey the equation 

(49) 



(J If -i 

+ {X + ib + 2a 2 )n = DR 1 cos(# - 

o9 



The solution for r\ which is periodic in 9 is of the form 



T\ = A cos( 



+ B sin( 



(50) 
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where 



.4 



B 



DR^X + ib + a 2 ) 
to 2 + (A + ib + 2a 2 ) 2 
DR x uj 



LO' 



(X + ib + 2a 



2\2 



(51) 
(52) 



We now consider the small angular perturbations. We can substitute p = S(r— a)[l/27rr+e/i] 
into the continuity equation and keep terms linear in e to get 



dh dh KRxcosje 
dt 86 2na 2 







Assuming the periodic solution is proportional to e^ x+lb ^ as above, one finds 

f 1 = C cos(# - 0) + D sm(6 - 0) 

with 

C 



DR^X + ib) 



D 



27ra 2 (u 2 + (X + ib) 2 ) 

DRxu 
2-Ka 2 {u 2 + (X + ib) 2 ) 



(53) 

(54) 

(55) 
(56) 



We can now rewrite the steady-state equations stemming from (31) and (32) in terms of 
the density to get 

pD + J 



-R 
pD 

b + u 
~p~D~ 



OO f2lT 



-oo JO 



rcos (9 — 4>)pr d6 dr g(u) du 



oo Jo 



oo p2ir 



rsin (9 — (f)pr d6 dr g(u) du, 



(57) 



where we have used that the order parameter for the solutions is chosen so that ^ = b 

ping term 

g(cu)du + 



Plugging in (46), (50), and (54), and keeping terms first order in e, 

X + ib + 2a 2 



2(pD + J) 



pD 

2{b + u ) 
pD 2 



X + ib 



(X + ibf + uj 2 



g(cu)du + 



(X + ib + 2a 2 ) 2 + oo 2 
{X + ib + 2a 2 ) 2 + u 2 



g(oo)duj (58) 
g(u)du (59) 



, (X + ib) 2 + u 2 

The bifurcation condition requires that A = 0. So the stability boundary is given by setting 
A = in the equation above. This gives (using usual relationships for principal values of 
integrals in the limit A = + ) 



2(pD + J) 
pD 
2(b + u ) 
pD 2 



ng(b) + 



ib + 2a 2 



P 



dbj 



oo(0 + 



UJ 



ib + 2a 2 ) 2 + u 2 



+ 



g{ui)dui 



to 



(0+ + ib + 2a 2 ) 2 + u 2 



g{uj)dio 



(60) 



(61) 



where P denotes the principal value. For the line D = Aq (i.e. a = + ) these equations 



reduce to (42) showing that the incoherence joins the corner of the death state. 
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IV. EXPRESSIONS FOR LORENTZIAN & RECTANGULAR DISTRIBUTIONS 



Lorentzian Distributions 



We now derive analytic expression for the boundary for the special case where g{oS) is a 
Lorentzian, 

9^) = \rrr772 ( 62 ) 
For a Lorentzian distribution, the equations are particularly simple because the Fourier 
transform is a simple exponential: 



g(p) 



-\p\r 



(63) 



We now plug this into the equations for the stability of amplitude death (25 ) and use the fact 



that these equations can be thought of as a convolution for b. A straightforward calculation 
then shows that the resulting equations for the stability boundary are identical to the case 



where g(oo) = 5(lo) } given by (23), except with D — X — > D — X + T, 

pD + J D - A + T 

pD 2 (D - A + T) 2 + b 2 
b + iojQ b 



(64) 



pD 2 ' (D - A + r) 2 + b 2 ' 
Thus T has the intriguing effect of decreasing the effective Ao, thereby pulling the individual 
oscillators closer to their supercritical Hopf bifurcation. 



B. Rectangular Distributions 

The integrals can also be performed for the case where g(u) is drawn from a rectangular 
distribution, 



g(u)=l/V if-r/2<w<r/2 
= otherwise 



In this case, the integrals in (24) and (25) can be performed and yield the equations 
a + A 



pD 2 



2T 



(arctan[(6 + T)/(a + B)} - arctan[(6 - T)/(a + B)]) 



b + u 
pD 2 



1 

2T 



log 



(b + T) 2 + (a + B) 2 
(b-T) 2 + (a + B) 2 



(65) 

(66) 
(67) 
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